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I. INTRODUCTION 



Optical lattices are gases of ultracold atoms trapped in periodic potentials created by 
periodically arranged intersecting standing waves of laser light. The interest in experimental 
and theoretical investigations of these artificial crystals is caused by the two following factors 

1) Neutral atoms in these optical lattices have several of attractive features that make 

n 

them interesting candidates for the realization of a quantum computer [2J. 

2) They may be used to simulate various lattice models of fundamental importance in 
condensed matter physics. Since they permit studying in a controlled way solid-state physics, 
in which one can fine-tune the interaction strength for various geometries of the lattices. In 
particular, it is possible to control the Hamiltonian parameters and study various regimes 
of system parameter. 

The lattice of bosons with short - range repulsive pair interaction trapped in an optical 
lattice may be described by a Hamiltonian of Bose-Hubbard type: 

tj n s N a 

H=-jj2btk+-Y: + E( £ i - a06&, (i) 

(iJ> > 1 

where b\ and b\ are the bosonic creation and annihilation operators on the site i; the sum 
over (i, j) includes only pairs of nearest neighbors; J is the hopping amplitude, which is 
responsible for the tunneling of an atom from one site to another neighboring site; U is the 
on site repulsion energy, and N s the number of sites. 

At zero temperature with an integer filling factor v = N/N s , where N is the total number 
of atoms, a system of bosons described by the Hamiltonian ([T]) could be on superfluid (SF) 
or in Mott insulator (MI) phase. Clearly the quantum phase transition (QPT) between these 
two phases is allocated by the dimensionless interaction strength parameter u = U/J. For 
small u, the hopping term dominates the system, so that it prefers to be in the SF phase. 
For large u 3> 1, on the other hands, the system exhibits a MI phase. 

A critical interaction strength w crit = 29.34 was found for d = 3 by Monte Carlo calcula- 
tions J3J at a filling factor v — 1, and agrees well with the experimental data (4]. 

To make for easier reading, we summarize some specific features of these two phases. 
The SF phase is characterized by a long-range correlation, a continuous (gapless) excitation 
spectrum and a finite compressibility. Since there exists a condensate with a finite number 
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of particles, no, the gauge symmetry is spontaneously broken in accordance with theorems 
by Bogoliubov and Ginibre. In contrast, in the MI phase, there is no long-range correlation 
or breaking of gauge symmetry. The excitation spectrum has a gap and the system is 
incompressible, since there is a fixed number of atoms per-site. The mobility of atoms is 
completely different in the two phases. In the SF phase they can easily move from one site 
to another site by tunneling, whereas in the MI phase, they are localized. 

Finite-temperature phases of optical lattices have been studied by Monte Carlo calcula- 
tions as well as experimentally for d = 3. As expected, the system behaves as a normal fluid 
(NF) at T > T c . A most intersting observation was made in Refs. 3j, |4|: In contrast to the 
system of dilute Bose gases, the critical temperature is downshifted at the transition to the 
MI phase. 

Most of present theoretical approaches can be divided into two classes which deal with 
the phase transitions SF— >MI and MI— >SF. The latter are based on the Ginzburg - Landau 
theory as described for instance in Ref . Q . They are well suited for analyses of time-of- flight 
pictures of the condensates and the associated visibility at zero and finite temperatures. Al- 
though such approximation works well in MI phase, predicting Mott lobes with quantum 
phase boundaries, they can hardly be extended to the SF regime in practice. In the former 
class (SF— j-MI) one uses a perturbative scheme [6| within a decoupling (or single site) ap- 
proximation due to Gutzwiller, or within a mean field approximation. Gutzwiller variational 
approach which was first proposed for a fermion system [7[ , and further developed for bosons 
in Refs. lii], has the following drawbacks [lol |: 

• The mean field Hamiltonian which features single boson terms does not conserve the 



total number of bosons 



Hi; 



• Tunneling of uncondensed atoms is neglected; 

• The critical value u cr - lt does not depend on the lattice dimension. 

The first application of a mean field approach was made in the Hartree - Fock-Popov 



(HFP) approximation to optical lattices by Stoof et. al. [12|. By studying the dependence 



of the condensate number n on u, i.e. n (U/J) they observed that no never reaches zero 
even in the unitary limit (u — > oo) implying that this approximation is unable to predict a 
QPT of SF — > MI. In contrast to this, the two-loop approximation by present authors in [3] 
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suggests the existence of such a QPT , but the critical value of u CI i t was found to be rather 
small: u cr ; t (two- loop) ~ 6 for d = 3. So, the question about the ability of approximations, 
based on mean field theory to describe phase diagrams of optical lattices adequately remains 
still open. 

A new approximation has been recently made by Fred Cooper et al. [bll . [ljj] in the 
treatment of dilute Bose gases. They found a simple way of fixing the degeneracy in the 
elimination of the interaction by auxiliary collective pair and density fields via standard 
Hubbard - Stratonovitch transformations. Although their approach gives no QPT for a ho- 
mogenous Bose gas at zero temperature, it predicts a desirable second order BEC transition 

t in the critical temperature T c that is 
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171 ] . One of the novel features of that 



at finite temperatures and exhibits a positive shif 
consistent with Monte Carlo an other calculations 
calculations is that for T < T c it predicts a novel type of superfluid phase that does not 
have a condensate [jjj]. Although such a phase has not been observed yet, it was justified 
by the existence of a nonzero anomalous density S, in the region T c < T < T*, where T* is 
the transition temperature to the normal phase . 

In the present work we shall formulate a similar two-collective quantum field theory for 
discrete systems such as optical lattices and ask the following questions 

• Does it predict a SF — > MI quantum phase transition? 

• Does it predict the suppression of T c at large u ? 

• Does it predict a new phase , mentioned above, for optical lattices either? 

Our results will be compared with those of another well-known mean-field approximation, 
the Hartree - Fock - Bogoliubov (HFB) approximation, which is widely used to describe BEC 
in homogeneous Bose gases and in triplons [lj], 0] in magnetic insulators, and will also be 
extended here to optical lattices. 

The paper is organized as follows. In Sections II and III we shall derive Collective 
Quantum Field Theory and HFB approaches for optical lattices, respectively. The results 
and discussions will be presented in Section IV, and the conclusions will be stated in Section 
V. 
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II. COLLECTIVE QUANTUM FIELD THEORY OF 3D BOSE-HUBBARD 
MODEL 



In the Wannier representation the Euclidian action, corresponding to the Bose - Hubbard 



13| 



Hamiltonian is given by 

A(r^) = / /3 rfr{^^(x i ,r)[9 r -/i]^(x i ,r)- J^^(x i ,r)^(x J ,r) 
Jo I i 



+ 77^*( Xi ' T )^*( Xi > T M X i> r M x i' r ) j> 



(2) 



where [i is the chemical potential and (3 = 1/T. The lattice points lie at the positions 21] 



Xi = i a, 



(3) 



where a is the lattice spacing, and 



i = (h,i 2 , ■ ■ .,id), 



(4) 



are integer-valued vectors. 

The partition function Z, and the grand thermodynamic potential Q, can be found as: 



Z = j Di)*Di)e~ A{r ^\ 
n = -TlnZ. 



(5) 
(6) 



The ground state expectation value of an operator O (ip*, ip) can be expressed as a func- 
tional integral: 



Z 



(7) 



With the help of a Hubbard - Stratonovich transformation, the interaction term in ([2]) 
can be eliminated by adding to the action in the exponent of ([5]) a dummy functional integral 
/ VAVA*e~ A ^ containing a pair field A with following action [22[: 



^pair = dT l2 j^/| A ( Xi ' r) " ^(x i ,r)V(x i ,r)| 2 | 



(8) 



A functional integral over A, A* produces an irrelevant multiplication of the partition func- 
tion Z by a trivial constant factor. 



It has been emphasized in 22| and the textbook 23| that this procedure is highly degen- 
erate. Instead of OH]) one could just as well have introduced a plasmon field by adding to 



the action in the exponent of (|SJ) the dummy functional integral containing a density field 
if with following action 



C/^(x i ,r)^(x i ,r)) 2 . (9) 



Diagrammatically, the degeneracy is caused by the fact that the sum of all collective 
field diagrams will always produce the same result if evaluated to all orders in perturba- 
tion theory. Each of these collective fields reproduces all effects of the interaction if it is 
integrated functionally. A difference appears, if the evaluation is restricted to a mean-field 
approximation. Then it depends on the dominance of certain dynamical effects which field 
is preferable. 

In principle, we can also add a combination of -4 pa ir and A p \, and still leave the physical 
properties of the system unchanged. For instance cosh 2 8A P \ — sin 2 6L4 pa ; r . Diagrammatically, 
however, the degeneracy cannot be easily verified since a calculation of the diagrams to all 
order is really impossible. It can only be done to some finite order, for instance in a loop 
expansion, so that the mathematical equivalence is initially of little use. 

One method to avoid the degeneracy and make the collective field approach unique has 
been pointed out a long time ago. It is based on an extension of the standard effective po- 
tential approach form from ty], whose expansion terms are the one-particle irreducible 
vertex functions of the theory. A unique version of collective fields can be introduced by 
going to a higher effective action A[^*,^f,A,A*,ip\. While the ordinary effective action 
w 4.[\J/* ) ijr] i s derived from a Legendre transformation of the generating functional of the the- 
ory W[r], 77*] in which additional source terms rjijj* + r]*ip have been added to the action, the 
higher effective action is obtained from the Legendre transformation of a generating func- 
tional W[t), rj*,j, K, K*] in which additional sources have been added to the action coupled 
to the density and the pair fields. The higher effective action will depend on the expecta- 
tions of the fields ip, ip*, p oc A oc iptp and A* oc At the end, it must merely be 
extremized, and no extra functional integrals 24| can cause any double- counting of Feynman 
diagrams. The expansion terms in the higher effective action are the two-particle irreducible 
vertex functions of the theory. 

Another method that also abandons the fluctuations of the collective fields in favor of 
a collective classical field has been dev elop ed in recent years from a generalization of a 
variational approach to path integrals 25| to all orders in perturbation theory. It was 
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extremely successful and has led to the most accurate theory of critical phenomena 26( so 



far, named Variational Perturbation Theory (VPT) (for a review paper see [271]). 

A third method which has recently been proposed and applied flQ uses the com- 
bination of both fully fluctuating collective fields implied by the above dummy action 
cosh 2 9A P \ — sin 2 9A pa i T for the particular value sinh^ = 1. This choice is preferable if 
we want the mean-field approximation to exhibit excitations that have no energy gap, to 
comply with the Nambu-Goldstone theorem. After a trivial change of the normalization of 
plasmon and pair fields in the total action A + cosh 2 A p \ — A va ,i r one arrives at 

A = A^*,tp] + A p [(p] + A A [A,A*), (10) 

with 

rP 

A^[ip*,ip] = / rfr^{^*(xi,r)[9 r -// + cosh0<p(xi,T)]^(xi,r) 
i 

3 sinh 0[A^*(x i; r)V(xi, r) + A^(x i; r)^(x i; r)]}— J f dr ]T ^*( Xi , r)^( Xj , r),(ll) 

^] = -{^E^ ^[A,M-f^E A(Xi ' r ^ (Xi ' r) - (12) 

At the level of for fully fluctuating fields if, A, A*, the parameter 9 is still arbitrary. The 
enforcement of the Nambu-Goldstone theorem will fix the parameter 9 to make sin 9 = 1. 
Now we consider separately two regions, with and without a condensed phase. 



A. Condensed phase 



In this phase, the U(l) gau ge s ymmetry is spontaneously broken. It can be studied after 



a Bogoliubov shift of the field 



io| 



(13) 



with 



(14) 



where the no = N /N is the condensate fraction. It is a constant in the absence of a 
magnetic trap. The fluctuating field ip(x, r) must satisfy the condition: 



f drX)fe,T) = 0. 

JO : 



(15) 



Substituting ( {TBI into ( |TTT) . and decomposing the quantum field ip('Xi,t) into its real and 
imaginary parts ^i(xi,t) and ^ 2 ( x i>0 as 

1 1 

^(xi,t) = -^=(^i(xi,t) + iip 2 (x h t)), ^*(xj,t) = -^=(^i(xi,t) - z^ 2 (xi,t)), (16) 

we may separate the action as follows: 

A = A + A 2 + Aa+A, p , (17) 

with 

A = -N s f3un (fi + Jz ) + im V / rfr [cosh 6(p(x u r) - § sinh 6> (A(x i; r) + A*fe, r))] ,(18) 

1 r _. - _ n 

A 2 = -J2 dr J2 \ie ab ip a {xi,T)d T ip b {xi,T) +ip a {xi,T)X a i! b (x U T)S ab (19) 

^ i J ° a,6=l,2 

~ o / ^EE^f^'^"^' 7 )- ( 20 ) 

where ^4a and ^ are given in f lT2|) . e a6 is an antisymmetric tensor with e 12 = —£21 = 
1, z Q = 2d = 6, and 

Xi = —fj, + ip(xi, t) cosh — \ sinh# (A*(^i, r) + A(x;, r)) , 

(21) 

X 2 = —\i + ip(x\, t) cosh. 9 + 1 sinh^ (A*(xi,r) + A(x- u t)) . 
For a homogenous, system the condensate is uniform and it is convenient to decompose the 



fluctuations into a Fourier series as |28l. |29| 

$ (x,, r) = —^= ' I ^a(q, tOn)e~ iuJnT exp [*x iPq ] (22) 

where u n = 2imT are Matsubara frequencies, and p q = {q±, q 2 , ■ ■ ■ , qd} 2n/N s a, with q\ 
running from 1 to N s — 1 are the discrete-valued momentum vectors in the Brillouin zone. 
The momentum sum is explicitly 

-, -, N S -1N S -1 N a -1 

jE'4^^-^' < 23 > 

JVs q ly s qi=0 g2=0 qd=0 

The prime on the symbol indicates that the p = -mode is omitted since it it contained in 
the subtracted ipQ. This will be useful to avoid possible infrared divergencies, especially for 
d < 3. 

In momentum space, the quadratic term A 2 reads 

A 2 = ^ ]T Va(q,^n)G~ 6 1 (q,o; n ;q / ,w m )'(/'b(q / ,a; m ), (24) 



q,q ,rn,n 
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with the propagator 
G(u n , q 



e(q)+X 2 -Jz uj n 
-u n £(q)+Xi-J2 



^ + ^( q ; 

where the bare dispersion e(q) and phonon dispersion £(q) are given by 

e(q) = 2j(d- 5>os(27rg a /Ag), 

^ a=l ' 



(25) 



£(q) = ^/(X 1 + £ (q)- Jzo)(X 2 + £(q)-Jz ). 

In the long- wavelength limit, e(q) behaves like 

Air 2 



(26) 
(27) 

(28) 



By comparison with the usual momentum-dependence of a free single-particle energy p 2 /2M 
we identify the particle mass M = l/2Ja 2 . 

Note that in coordinate space the Green function is defined by 



G ab (xi, r; xj, t) = G ab (xi - xj, r - r) = (^a(x i; r)V>&(x.j, r )) 
= E E e^^V^i-^K, q). 



(29) 



The thermodynamics of the system can be calculated from the partition function Z func- 
tional integral over all fields ipi, ip 2 , V 9 ; A and A* fields 

Z = J VTp 1 Vip 2 VipVAVA*e~ Ao ~ A2 ~ AA ~ A ^. (30) 

The first integrations by ifii and ^2 are Gaussian and may be evaluated easily by using 
well-known formula 



2 E f^a(x)Gj(x,y)i) b (y)dxdy-Jj 1 (x)i) 1 (x)dx-Jj 2 (x)'i()2(x)da 



a,b=l,2 



cxp 



a,6=l,2 



(31) 



The integrations over the fluctuating collective fields, however, cannot be performed exactly, 



since they are nontrivially contained in VDetG. As usual in these circumstances, we resort 
to the saddle-point approximation 22|, |30|. In the absence of a trap, we may assume the 
saddle point to lie at constant values of ip{x^r) and A(xi,r): 

<p(xi,r) = ipo, 
A(x;,r) = A*(x;,t) = A . 



(32) 
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Then the integrals over ip a become trivial and we may use the formula Det G = e Tr ln G in 
Eqs. fl30|) and ( l3~lj) to derive the following effective potential: 

with 

A = Asinh#, cp' = —fx + cosh#(/2o — J^o- (34) 
The spectrum of density fluctuations is now from ( 1271) : 

£ 2 (q) = (e(q) + - A)(e(q) + <p' + A). (35) 
The sum over p may be calculated in d = 3 by approximating f[2"3l as follows 

ir £ /(e(q)) -> / dq 1 dq 2 dq^{e <l ), (36) 

with the lattice dispersion: 

3 

e q = 2j]T [l-cosTrgJ. (37) 

Note that on lattices, the momentum integrals are always finite so that there is no need 
for renormalizing the coupling constant. This is in contrast to atomic gases. However, if we 
want to express the coupling constant in terms of the scattering length a s that is observable 
at low-energy in atomic gases, where the quadratic coupling constant g must be renormalized 
to a finite value qr by the addition of a diverging integral 1/qr = 1/g + / d 3 p/ (2iih) 3 e(p), 
the relation a s = Mgn/Anh 2 can only be employed only after a corresponding addition of a 
finite sum [see the remarks after Eq. (1921)]. 

Another remark concerns the frequency sum in (1331) . which is is initially divergent. In 
fact, to evaluate a frequency sum such as Y^=-oo l n (o 2 + w n) with u n = 2iinT, one must first 
differentiate it with respect to a, perform the summation over n, and integrate the result 
over a 31[. This procedure gives an additional divergent constant, which may be removed 



by an additive renormalization of the energy 32j. The subtraction can actually be justified 
by calculating the path integral as a product of individual integrals for each slice of a sliced 
time introduced in the very beginning by Feynman [31] . 

In the thermodynamic potential Q, one subtracts from Q the one for the "ideal" value 

Sl(U = T = 0) = \ $>(q) -fi- Jz ) = \ 5>(q) + </), (38) 

9 Q 
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and deals only with the subtrated expression 



Q ren = Q{U, T) - n(U = 0, T = 0) = - £(£(q) - e(q) - tp 1 ) + N s vn ( V ' - A) 

+ NA \ - ^ + " + J *> 8 + T £ ln(l - e-^W), 
2t/sinh 2 # 2f/cosh 2 # V 

where we have performed summation by Matsubara frequency by using formula 

oo 

^2 ln(u; 2 + a 2 ) = a/3 + 2 ln(l - e _/3a ) + divergent constant. 

n=-l 



(39) 



(40) 



For brevity, we shall suppress writing down the subtraction in f2 rcn . 

In equilibrium, the thermodynamic potential reaches a minimum with respect to para 
meters no, <p' and A. Thus we minimize fl with respect to n 

on 



dn 



N s u(<p' - A) = 0, 



and get 



ip' = A. 

Inserting this into ( 135|) leads to the well-known Bogoliubov phonon dispersion 



£(q) 



(qW<-(q)+2A, 



(41) 



(42) 



(43) 



which is linear in q for small momentum, thus respecting the Nambu-Goldstone theorem. 
Minimizing thermodynamic potential fl with respect to A gives the equation: 



A = U s'mh 2 



A ^ C q 



(44) 



where c q is short for 



c q = \ + //»(£(q)) = \ ™th (fiE (q)/2), /^(w) = l/(e*" - 1). 



(45) 



Minimizing with respect to ip', thereby taking into account the relation d£(q)/dip' 
(e(q) + v?') / f (q) , gives the following equation: 

N s (p' + + Jz Q ) 



N s vn + 



(e(q) + yQcq 1 
£(q) 2 



C/cosh 2 6» 



0. 



(46) 



This will serve to determine of uncondensed fraction n, 
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B. Normal and anomalous densities 



According to the general rules of statistical mechanics, the total number of particles N 
is conjugate the chemical potential: 



N = - ^ 



Applying this to ( )39l) gives 



N 



N s {ip' + n + Jz ) 
[/cosh 2 (9 



Using ( 1481) in ( 1461) . we obtain 

N = N s un + 



(e(q) + (p')c c 



1 



N + N u . 



(47) 



(48) 



(49) 



f(q) 2. 

Here N is a total number of condensed atoms, and no = Nq/N s v is the condensate fraction. 
The uncondensed atoms have a fraction 



Nu 

N 



vN, 



AT 



(g(q) + <p')cc 



(50) 



It satisfies the trivial relation no + n u = 1. 

Note that, the term — | in the square bracket of (1501) is due to the renormalization 
procedure ( 1391) . and guarantees that at T = all particles of the ideal gas (which has U = 
and A = 0) are condensed, so that n-JJJ = 0, T = 0) = 0. 

If the U(l) gauge symmetry is broken, a Bose system is characterized not only by the 
expectation values of the fluctuating part of the ^-field with the normal density n u = (tp*tfj), 
but also with anomalous density, defined by 



8(xi,T,Xj,T') = (^(x i ,T)^(x i ,T / )). 



(51) 



Clearly, for homogenous system in the equilibrium, in particular, for periodic optical lattices 
without magnetic trap, 5 does not depend on coordinates, i.e. 5(xi, t, Xj, t') =const as was 
emphasized in [33]. Omission of the anomalous averages makes all calculations not self- 
consistent, the dynamics non- conserving, the thermodynamics incorrect. It ruins the order 
of the phase transition and renders the system unstable. It was also shown in j^J that a 
5 = type of mean-field approach referred in the literatures as Hartree-Fock-Popov (HFP) 



approximations 



19| leads to a discontinuity in the magnetization curve of antiferromagnetic 



material with triplon BEC. Thus we must always allow for 5 ^ 0. 
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Let us calculate this expectation value from the formula 

S = -(^(xi,r)^(a;i,r)) = ^[(^i(xi,r)^i(xi,r)) - (^ 2 (xi,r)^ 2 (a;i,r))] 
v Iv 

= ^[Gii(0)-G? a (0)]. (52) 
In momentum space, the propagator can be rewritten as 

GotK , q) = __^ ^j, (53 ) 

where we used equations fl2~T|) . fl25|) . and f|43l . Using in f )52|) the equations fl29|) and fl53l) . 
one obtains 

1 2A A v c q A v 1 /l 1 \ 

2 ^n q^n+^(q) " N * V ^ q ^ ^ ~ 1 ' 

In terms of 5, the A-equation may be rewritten in the following compact form 

A = Uv(n + 5) sinh 2 (55) 



with rio = 1 — and n u given by ( 1501) . 

It is well known that the Goldstone theorem for a dilute Bose gas with a spontaneous 
broken symmetry is equivalent to the celebrated Hugenholtz-Pines theorem 34J, according 
to which self-energy S c i and the anomalous self-energy A c i satisfy 

£ci - A ci = (56) 
In the Appendix A we shall show that a similar equation holds for optical lattices: 

E c i - A c i = H + Jzo, (57) 

with S c i = ifo cosh 8, A c i = A. 

The only parameter, that so far remains free in the initial action (fTTj) . is 8. It may be 
chosen such that the quasiparticle energy £(q) reduces, in the one-loop approximation [13], 
to the gapless Bogoliubov dispersion 



£ (q)oneioop = V e (q) V 6 ^) + 2Uv - ( 58 ) 
Indeed, in this approximation we get from ( )55l) A m Uusmh. 2 8, and from ( )43|) £(q) ~ 



e(q)ye(q) + 2£7^sinh 8. And this is the place where we fic the 8 to satisfy 

sinh 2 9 = 1, cosh 2 6 = 2, (59) 
13 



as was announced earlier. 

Summarizing this section, we present the full expression for Q: 



n = -J2[£(q)-e( q )-A} + 



N S A 2 N s (A + n + J 



^2 



AU 



with 



H = 2ull - A - Jz . 



+ T^ln(l-e"^ (q) ), (60) 



(61) 



The last equation follows from ( )48l) . The self energy A in ( l60j) and ( l6~il is defined through 
the following set of nonlinear algebraic equations: 



A = Uu(n + S), n = l-n u , 



E 



c q (e(q) + A) 1 
£(q) 2 



(62) 



WV s ^£(q)' 

where c q is given in ( T45i) and £7, J, i/, T are input parameters. 



C. Symmetric phase 

When no = 0, the Hamiltonian (CQ) is symmetric under the transformation ip — > e^ip and 
equation (T4TT) makes no sense. Then 7^ A, and the energy spectrum has a gap with the 
dispersion 



£(q) = ^(e(q) + ^ - A)(e(q) + iff + A). (63) 
The main equations in this regime with T > T c are 



A = TJvL 5 



A 



WV S S f (q) ' 17 jv fl E 



(e(q) + ^)cq 1 



(64) 



5(q) 2 

The set of equations ( 1641) with the energy spectrum ( l63i) may have a solution A 7^ 0, <p > A 
leading to an exotic state with no condensate but with a finite anomalous density: n = 



0, 5 7^ 0. It was shown in Ref. [18] that this phase has a nonzero SF transition. The upper 
boundary of such a state was denoted by T*, and was determined by solving the equations 
( 164]) with A = 0,tp' > 0. Thus it was predicted that dilute atomic gases posses a state with 
T c < T* that a superfluid but has no condensate. However, experimentally, the states with 
T > T* are normal with 5 = and no = 0. In Sect. IV we shall investigate the possible 
existence of such a state for optical lattices, with a negative outcome. 
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III. VARIATIONAL PERTURBATION THEORY IN OPTICAL LATTICES 



Let us compare our result with those of Variational Perturbation Theory 26j] . To lowest 



order, this is equivalent to the HFB approximation used in the operator formalism |35j. To 
do this, let us formulate the HFB approximation for optical lattices in the functional integral 
framework. 

Starting point is again (Tj[|) in which we perform the Bogoliubov shift ffl3|) and separate it 
as follows 

A = A {0) + A {1) + A {2 ) + A {3 ) + .4(4) , ( 65 ) 

where 

.4(0) = f5N s vno[^vn Q - /x - Jz ], 

4.(1) = v / ^o[-A* - Jzo + Uun ] J dr^(^( Xi ,r) -f V>*(x ; , r)), 

fP ~ ~ 1 U 

•4(2) = j Q dT |z]V'*(xi,r)[5 r - /x]^(x i5 r) + 



Ei [^ 2 (xi, r) + #*( Xi , r)^(xi, r) + ^(x ; , r$*(xi, r)l - jV^*( Xi , T ) \ » 



(66) 



_ 

-4(3) = C/v / ^ y dTy^[ 1 //( Xi ,r)'0 2 (xi,r) +^*(xi,r)^*(xi,r)^(xi,r)], 



^(4) = Y y o ^E^*( x i' r )^( x i' T )l 2 - 

After this we add and subtract following terms 

As) = dT J2 {s c i^*(x i ,r)^(x i ,r) + ^A d [f ( Xi ,r)f ( Xi ,r) + ^(x i ,r)Vi(x i ,r)] J , (67) 

with variational parameters S c i and A c i. The subscripts cl emphasize that these are vari- 
ational parameters which, in contrast to the earlier fields (f and A, are not meant to be 
functionally integrated. 

Using again real and imaginary parts of the complex fields if), ip* as in ffl6|) . we rewrite 
A as 

A = A + A {rcc + A int , (68) 
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where 

1 [ p 

Aree = ^ ^ dr ]C X! Va (xi , r) [ie afe 5 T + y" a (5 a b]^fe (x ; , r) , 



1 /■/* 

2 Jo 

|( 2 ) J_ /I© J_ /l( 4 ) 



i a,fe=l,2 



•/Hint — «^int ^ "^int ^ «^int ' 

i r@ 

Ant = 2 i {^i(xi,r)[3[/^n - £ c l - A d ] + ^|(x i5 r)[f7z/n - £ d + A d ]|, 



i 

1 /"^ 

42 = -*7v^W / V[Vi(Xi,r) + Vi(x i ,r)^(x i ,r)], 

2 Jo j 

= \uj^drYX^,r) + ^ 2 (x i ,r)] 2 , 



where 

Yi = -/i - Jz + S c i + A c i, 

(70) 

y 2 = -m - Jz + S c i - A d . 
The free part of the action, Aree in Eq. ( |69|) . gives rise to the propagators to be used 
in perturbation expansions. In the momentum representation of the fields Eq. f l22|) . the 
propagator is given by 

/ 



G{co n , q) 



(71) 



with £ 2 (q) = (e q + Yi)(e <i + ^2)- To lowest order, one obtains 

= -T In Z = -T In Z - T In Z frcc + T(Ant) , (72) 

where Z = e"-\ Z hcc = j V^V^e^- = 1/VDet G" 1 , (Ant) = / V^V^A^e'^ . 
Now we evaluate 

(^ 2 (x i ,r))=G aa (0) = ^, (^(x i5 r)) = ^f, (Vftxj, r)^(x,, r)) = (73) 

with 

and we find the following thermodynamic potential: 

Q = iV s z/no ^-/i - J^o + "Tj-^o) + T^X^fa) ~ £ ( ( l) + /-* + ^0] 
+T$>(1 - e-^W) + ^[Sa 2 + 3a 2 + 2^] 

+-<7i(3[/i/n - F - J^o - A*) + ^(Uvrio -Y 2 - Jz - fji), (75) 
where we have again subtracted f2(T = 0, U — 0). 



1(3 



The parameters S c i and A c i are now determined variationally by requiring that they 
minimize the thermodynamic potential, i.e., we require dQ/dY> c \ = and dQ/dA c \ = 
or equivalently 



wr ' wr°- (76) 



These equations yield 



Y 1 = 2>Uvn - fi - Jzq + 2Wj(3cri + ^2), 

Y 2 = U vtiq - [i - Jz + 2^7(01 + 3cJ 2)- 
The gaplessness of the energy spectrum is now imposed by hand. In fact, by requiring the 
relation ( 157)) . we get from (170)) Y 2 — which leads to the dispersion 



f(q) = v / e (q) v / e (q)+2A, (78) 



where A = Yi/2. This leads to the equations 

A = Uvn + r^-Oi - 02), 

2Ns u (79) 
H + Jz = Uvn + ^-(ci + 3(72)- 

Here, we draw the reader's attention to the inself-consistency of the HFB approximation as 
far as the chemical potential is concerned. In fact, the stationary condition dQ/dn = 
with Q given by ( 1751) leads to the following equation for /1: 

li + Jz = Uvtiq + 7^f( 3cr i + 02), (80) 



which contradicts fi of Eq. (179)) . 



37l | proposed to 



To make the theory self-consistent, Yukalov and one of the authors 
introduce two chemical potentials: namely, yUo, which corresponds to the Eq. (180)) . and /ii 
corresponding to Eq. ( 1791 . Being responsible for subsystem of condensed and uncondensed 
particles they, naturally, coincide in the normal phase, when Yi = Y 2 = 0. In the present 
work, however, we follow the standard procedure of identifying /x in (179)) as a chemical 
potential from which we determine the particle densities by differentiation of Q. 



A. The fractions n u and 5 in VPT 



Applying the well-known relation N = —dQ/d/j, to Q in (175)) gives 

■(e(q) + A)c q 1" 

q 



N = N s un + 



£(q) 



n + n u , 



(81) 
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- !bi _L_ v 



(82) 



and hence 

- (e(q) + A)c q 1 
£(q) 2 

with the £(q) is the Bogoliubovs dispersion given in (ITS]) . 
For the anomalous density 5 we obtain 

5 = ife*,r» = ^L- [Gll(0 )- G22 (0)] = <°0> 

= -wv:?fw- (83) 

where we used Eqs. (ITIjl and (174|) . 

Using now ( |83l) in ( 1791) gives the equation: 

A = Uv(n + 5), (84) 

which is the same as the one in before (155]) with fl59l) . The only difference between these two 
approximations is in the sign of anomalous density, which is 5 > in the collective quantum 
field theory and 5 < in HFB. 

Summarizing we collect here the main equations in both approximations: 

A = Uu{n + 5), n = 1 — n u , (85) 

* = £^E^ 5(q) = v^(q)^(q) + 2A . c ^ + ^nrr' ( 86 ) 

p = 2Z7z/ - A - Jz , (87) 
f -1, HFB , 

e = (ss) 

I +1, Two Collective Quantum Fields. 

Note that similar relations hold for atomic gases. A difference occurs for the T > T c 
phase. There one must use replacements listed in Appendix B. In fact, in the normal phase, 
n = 0, HFB theory gives 

Since the right-hand-sid of this equation is positive, while the left-hand- side is negative, at 
least for optical lattices, Eq. (18"9|) has as only solution A = [see Eq. (155|) ]. This means 
that in the normal phase no = and 5 = simultaneously. Therefore HFB theory does not 
predict a superfluid phase without a condensate, this being in contrast to the two-collective 



quantum field result of Cooper et.al. in Ref. 18J at the mean-field level. 
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From above discussions it is easy to understand that VPT gives no shift in T c due to 
interaction. In fact, when T — > T c , the condensed fraction no — > 0, and hence A — > 0. The 
expression for n u , will coincide with that for the ideal gas, i.e., Eq. becomes 

V = ? - i = iv, ^ e^)/ 7 ? - 1 ' (90) 

which means that T c = T c ° for HFB and, hence, AT C = T c — T c ° = 0. 

IV. RESULTS AND DISCUSSION 

A. Quantum phase transition in two-Collective Quantum Field Theory and VPT 

First we discus the existence of QPT in optical lattices for two collective quantum fields 



at the mean field level and for the HFB approximation. It has been shown that 



atomic Bose gases Collective Quantum Field approximation does not predict QPT 15] while 



or dilute 



HFB does [33|. Below we show that in the case of 3d optical lattices the situation is vice- 
versa. This can be understood in the following way. Lets rewrite main equation at T = 
as: 

no(A) = A_$(A). (91) 

It is clear that for interacting system, U ^ and A ^ 0. Since in the collective quantum 
field theory 5(A) > 0, the Eq. (19T|) may have solution n (A) = with A ^ (see Table 
I). However, in HFB approximation 5(A) < and no (A) in (|9T|) may have only solution as 
n > for A ^ 0. Note that in the case of dilute atomic gases 



5(A) = —8py/j/ir < two-Collective Quantum Field Theory 
5(A) = +8p v / t7^ > HFB, 

with the dimensionless gas parameter 7 = a^p that characterized the interaction strength 
of the gas after renormalization. It is formed from the s-wave scattering length a s and the 
particle density p. This sign change is responsible for the dilute atomic gases has a QPT in 
the HFB approximation, but not in the two-collective quantum field theory at the mean-field 
level. 
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TABLE I: Critical parameter u c = (U/J) c of QPT in the two collective quantum field approach. 
The critical temperatures of ideal optical lattices in d = 3 are listed (in units J). The third row 
presents approximated values of t° c [see Eq. 



V 


1 


2 


3 


4 


5 


u c = (U/J) c 


28.04 


47.70 


67.15 


82.52 


105.91 


4 = Til J 


5.6 


9.69 


13.70 


17.70 


21.67 


t® in small q approximation 


5.06 


10.07 


15.2 


20.25 


25.32 



B. Critical temperature T® for ideal cases 

Before we study the shift of T c , let us estimate the critical temperature T c ° for the free 
optical lattice with U = 0. Assuming A = in Eq. (IMj) . we obtain the well-known formula 

(93) 



dq 1 dq 2 dq 3 -^Q- 



Introducing dimensionless parameters t 
may rewrite f )93|) as 



e- 
T C /J, 



1 



£ q = ej2J = ELil 1 - cosvrg a ), we 



i 1 
dqidq 2 dq 3 —- W) - 



(94) 



which can be considered as a nonlinear equation for t° at a given filling factor v. Our 
numerical estimations for t° c are given in Table [U It is seen that for v — 1, T° = 5.6 J, which 
in consistent with other estimates given in the references j^, [h]] . 

Note that T c ° can be approximated as T c ° = 5.6 Ju - 825 in the range v G (1,5) including also 
non integer values. In the third row of Table 1 an approximated values of t° are presented. 
This approximation, say, spherical approximation at for small momentum is obtained by 
following replacements in (|94|) : 

■i 7r rid _ 



dqidq 2 dq 3 f(q) - / q'dqf(q), 
2 Jo 

2 



2 {(fq /n _ 



7^0 

c 

e(q) 



where the Debye momentum grj defined by the equation: 

it r qd 



l 



dqidq 2 dq 3 = - 
o 2 Jo 



q 2 dq, 



(95) 



(96) 



equals to ga = (6/vr) 1 / 3 ps 1.24 for d = 3 which gives T c °/J = 2z/7r(7r/6) 1 / 3 . It is seen that 
this approximation works with roughly 10% accuracy for v < 3. 
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C. The shift in T c caused by the interaction 



We are now prepared to estimate the shift AT C /T C ° = (T c — T c °)/T c ° analytically. The 
integral is dominated by small momenta. The main equations at T — > T c for n = 0, n u = 1 
are 

A = UA f'dq^dq M 8 .^ , (97) 
Jo S(q) 

1 = - f 1 d qi dq 2 dq 3 ^±p- f B (E(q)), (98) 
v Jo t{q) 



with£(q) = ^Ve q + 2A, / B (f(q)) = l/(e^-l), & = 1/T C . 

Note that in fl97|) we may assume A ^ and divide both sides of f l97|) by A. The critical 
temperature of ideal gas T c is the solution of Eq. ( 198!) with A = 0, i.e., 

1= W'W» (99) 

f jo e £c > / - , c — 1 
Now we introduce dimensionless variables: 

A = u 2 k 2 T°, T c = T>, T c ° = Jt°, e q = 2Je q , £(q) = 2J£(q), (100) 



with e q = ~~ cos7rq a ), £(q) = yt q yfc q + u 2 kH®, AT c /T c ° = a — 1 and t° c are given 

in the third row of Table 1. 

The scaled equations can be rewritten as follows: 

= j _ w Z" 1 fB(£(q))dqidq 2 dq 3 

2 Jo ' £(q) ' 1 J 

= 1 - £ f dgi^ 2 dg 3 gq+ ^ /2 /B(g(q)), (102) 

2 jo £(a) 



with / B (£(q)) = l/( e 2 ^(q)Mc _ i). 



Bearing in mind (|99|) . we may rewrite HI 02 j) as 

1 A^^S l — ^ 1 = 

^ ^ M \ e 2£(q)/tg-l 5(q)( e 2£(q)/or*g _ 1) J V 7 

The nonlinear equations fllOip and fl 1 U3j) should be solved with respect to k and a with 
given numbers u = U/J and t° c . To do this we make replacements fl95|) . Then Eqs. fllOip 
and (11031) can be rewritten as 

uat® r £ o de 



0, (104) 

1 1 I ..2. ,2,0 = °» ( 105 ) 



4V2tt 2 V v^(e + M 2 K 2 t°) 
cfe f a(e + u 2 kH° cI 
o yfe \ e + u 2 K 2 t Q c 
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where e D = 7i 2 ql/2 = (n 2 /2)(6/n) 2 / 3 

The integrals in (11041) and (11051) are easily done and yield 

= v / 2(6tt 2 ) 1/3 (1 - a) + uarey^ arctan (106) 
= 4vr 2 /t- v^ay^arctan^, (107) 

where 9 = v^tt 2 ) 1 / 3 / '(2«u^/*°). Excluding a from ffTOTl) and inserting it to ( 1106j) gives 

a = —j= -, (108) 

't9 arctan # 



= 4kvt 8/3 6 1/3 - ^2^ c [(6ti 2 ) 1/3 + 2uk 2 tt 2 } arctan^. (109) 
Now we consider separately two regimes: 

a) Weak interacting regime. Expanding f 1 1 8 p and (11091) in linear order by u we get 

4tt 2 k\/2 8k 2 m/6\ 2/3 
a = + - , (110) 



111) 



/ t o 3 W 



8tt 



Now inserting k into ( jl 1 0[) we finally obtain 

ut° / 6 \ 2 / 3 



and hence 



« = 1 + ^(^4) + 0(« 2 ), (112) 



AT . ui9 / 6 \ 2 / 3 



c 



which means that for small coupling constant, i.e. (U / J) < 1, the shift is positive. 

b) Strong interacting regime. In this region, A/u 2 i.e. k is small, so we may use a 
linear approximation in k in Eqs. (1 1 8 j) . (11091) 

^ 2 kV2 

a = p= — , 114 

= v/ ^ ( ^ 5)1/3 - n[2ut° c + 4(6vr 8 ) 1 / 3 ]. (115) 



This leads to following equation 



2^ 3 6 1 / 3 % 
01 ntp + 2^/36 1 /3 TP' 1 j 
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from which one may conclude that T c decreases with increasing u, i.e. 



AT C 



ut 



a 



< 0. 



(117) 



T c ° ut° c + 27r 8 / 3 6V3 

Thus, our analytical estimate shows that the critical temperature T c as a function of the 
coupling constant U, i.e. the function T(u) first increases and then decreases with increasing 
u for optical lattices. This is in good agreement with experimental measurements (4). Note 
that, this behavior is in contrast to the case of dilute atomic gases, where the shift is always 
positive in the whole range of repulsive coupling constant. 

In Figure 1 we present T c (in unit of J) vs. u for v — 1. The solid line correspond to 
the exact numerical calculation, i.e., the numerical solutions of Eqs. (|97|) . (|98|) . The experi- 
mental points (circles) are taken from 4|, solid triangles are from Monte-Carlo calculations 
taken from ref. 



10-, 




A Capogrosso-Sansone et. al. 
o Experiment 
O QMC 



I — 
10 



- 1 — 

15 



— I — 

20 



— I — 

25 
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FIG. 1: Behavior of T c as a function of U/J in the saddle-point approximation to the two-collective 
field theory for v = 1. The circles show experimental values given in [J], solid triangles are from 

Q 

Monte-Carlo calculations of Ref. [3:]. Note the initial rise that was found also in atomic gases in 
Ref. Lj 



The suppression of T c at large coupling constant is found for v > 1 as it is also seen in 
Fig. 2. 

In Fig. 3 we present the critical values of the self energy A c = A(T = T c ) in units J vs. 
(U/J). Observe that when J is fixed, A c increases with increasing u and v. On the other 
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hand we observed that A c in units T c °(z/), i.e. A C (T = T c )/T c °(z/) vs. u is almost independent 
on v, e.g., (A c /T°)\ v=1 = 10.425 and (A c /T°)\ u=4 = 10.447 at u = 21.0. 

Now we consider the behavior of A for T > T c . It was suggested by Cooper et.al. [3] that 
in the temperature range T e (T c , T*) there exists a ?7(l)-symmetric phase with no = but 
8^0. This would imply the existence of a superfluid state without a condensate. However, 
by solving (1641) for A and </?', we could not find, for optical lattices, any solution with A ^ 
0, <p' 7^ 0. Instead, the equations for T > T c , have a solution with A = 0, (p' — 2Uv— Jzq — u. 
In this normal state with 6 = 0, the filling factor that characterizes the particle density, is 
determined by the well-known equation 

1 1 f 1 1 

with the bare dispersion £ q = 2J^a=i(l — cos7rg a ). 

The chemical potential of interacting bosons in T > T c may be evaluated self consistently 
from Eq. (I118p with input parameters is, J, U , and T, or given by an external field (pumping) 



as in the case of triplons [19|, |20 |. 



V. CONCLUSION 

In this paper we have developed a Collective Quantum Field Theory and a HFB approxi- 
mations for d = 3 optical lattices at very low temperatures. We have shown that, in contrast 
to the situation in dilute atomic gases, a two-Collective Quantum Field treatment in the 
saddle point approximation predicts a Quantum Phase transition, that is missed in a HFB 
approximation. The shift of T c due to the point interaction is zero for HFB, while it has 
a nontrivial dependence on the coupling strength (U / J) in the Collective Quantum Field 
treatment. The latter agrees with experimental observations as well as with Mont-Carlo 
calculations. 

We have found no exotic superfluid state with finite anomalous density but zero conden- 



sate. Therefore, the temperatures T* and T c introduced by Cooper et.al. [18[ coincide. The 
system is in superfluid state for < T < T c , and in normal state for T > T c . Application of 
present theory to d = 1 and d = 2 lattices is in progress. 
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Appendix A 



Below we derive Hugenholtz-Pines theorem 

£ cl -A cl = / u + Jz , (119) 

of Sect. II for optical lattices. The normal S c i,and anomalous A c i self-energies in f II 1 9 j) corre- 
spond to the normal G Q (r,r') = (T r \l/(r)\l/ + (r')) and anomalous G an (r, r') = — (T T \I/(r)\l/(r / )) 
Green functions respectively. In the field representation (1161) we have: 

£ cl = i[n n + n 22 ], (120) 
A cl = i[n 22 -n u ], (i2i) 



where U a b are defined by Dyson-Beliaev equations 

(G^U - (G^U = n ab , (122) 
and the Green function Go corresponds to the noninteracting situation 

n - 1( v [ e(q)-/x-Jzb -u rl \ 

G (u n , q) = . (123) 

y u n e(q) - n - Jz J 

The interacting Green function G~ l is defined in Eq. (125]) . Using ([HJ, ([34]), ([25]), ffT23]) in 
ffT22l) gives: 

IIii = + = cosh^(/?o — A, 

n 22 = X 2 + /i = cosh #(£o + A, (124) 

n X2 = n 2 i = o. 

Inserting f lT24"]) into f lT20]) and ffT2T]) one derives 

S c i = <A) cosh 6*, 

(125) 

A d = A. 

and hence 

E c i — A c i = ipo cosh 6 — A = 9?' + /i + J^o — A, (126) 

where we have used Eq. (]34"]) . As it has been shown in Sect. II, in the condensed phase 
tp' — A and Eq. (ll26p becomes equivalent to the Hugenholtz-Pines theorem, i.e. to Eq. (II 19j) . 
The relation f 1 1 1 9 j) in HFB approximation can be proved in a similar way. 
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Appendix B 



Here we present formal equivalence between Bose-Hubbard Hamiltonian ([I]) in Wannier 
representation and standard Hamiltonian for homogeneous dilute atomic gases 

V 2 



H 



dr\pt(i 



2m 



*( r ) + § / rfr[*t( r )«r( r )]5 



(127) 



where g is the constant of contact interatomic interaction. Using the replacements listed in 
Table II we obtain for Q and the extremality equations in dilute atomic gases versus optical 
lattices the relevant quantified as derived in Sections II and III. Of course, an appropriate 
renormalization procedure is implied in dilute atomic gases. 



TABLE II: Formal similarity between Hamiltonian (TTJ and (|127p 



Quantity 


Atomic gases 


3D optical lattices 


Comment 


Volume 
Density 

Bare dispersion 

Chemical potential 
Momentum summation 

Normalization of densities 


V 

p = N/V 
e(q) = q 2 /2m 

M 

V/(e(q))=49r / q 2 dqf(e(q)) 

P0+ Pl= P 


N s 
v = N/N s 

3 

e(q) = 2 J S ^(1 — cos7rg a ) 
a=l 

p. + Jzo 

y^/(e(<l)) = / dq 1 dq 2 dq3f(e(q)) 
"0 + "1 = 1 


Af,- number of states 
filling factor 

no additional magnetic trap 
^ = -(f) 
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